function rie = computeRie(i,elements,model)
e = model.getElement(i-1);
f = e.getFace();
dim = 2;
nn = f.countModes();
ng = e.getNumGaussPoint();
rie = zeros(dim*nn,1);
for j=1:ng
    xi = e.getGaussPoint(j-1);
    w = e.getGaussWeight(j-1);
    jac = e.getFace().jacXAt(xi);
    detJac = det(jac);
    B = zeros(3,dim*nn);
    B = e.computeB(f,B,xi);
    sigma = elements(i).stress(:,j);
    rie = rie + w*B'*sigma*abs(detJac);
end
H = elements(i).H;
Gamma = elements(i).Gamma;
h = elements(i).h;
rie = rie - Gamma'*(H\h);
end
